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Abstract 

A  time  domain  method  is  presented  for  analyzing  simultaneous  measurements  of  pressure  and 
the  horizontal  components  of  velocity  obtained  beneath  irregular  multidirectional  wave  fields.  This 
new  method  differs  from  the  usual  linear  directional  analyses  applied  to  PUV  data  in  two 
important  aspects.  First,  the  essential  nonlinearity  of  the  measured  waves  is  not  sacrificed  to 
achieve  a  solution.  Therefore,  predictions  of  sea  surface  elevation  and  directional  kinematics 
throughout  the  water  column  accurately  portray  the  actual  nonlinear  character  of  the  waves. 
Second,  the  analysis  method  is  ‘local’  in  that  it  can  be  applied  to  segments  of  PUV  time  series 
much  shorter  than  an  individual  wave.  The  viability  of  the  locally  nonlinear  methodology 
developed  in  this  paper  is  proven  by  demonstrating  agreement  with  higher-order  theoretical  steady 
waves.  Predictions  of  sea  surface  elevation  and  wave  kinematics  are  also  made  using  actual 
measurements  from  PUV  instruments  at  two  ocean  sites  off  the  west  coast  of  the  United  States. 
©  1999  Elsevier  Science  B.V.  All  rights  reserved. 
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1.  Introduction 


The  PUV  gauge  has  been  a  popular  field  instrument  for  sensing  the  directional 
characteristics  of  waves.  It  groups  together  a  pressure  transducer  and  a  directional 
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current  meter  in  an  assemblage  such  that  both  instruments  have  the  same  horizontal 
position  and  adjacent  but  wholly  submerged  vertical  positions.  In  shallow  waters,  this 
assemblage  would  commonly  be  deployed  near  the  bed.  In  deeper  waters,  deployment 
may  be  much  higher  in  the  water  column,  because  of  the  considerable  attenuation  in 
signal  magnitudes  with  depth.  Because  these  gauges  are  not  surface-piercing,  they 
escape  the  most  extreme  environmental  conditions  and  have  the  potential  to  optimize 
data  recovery. 

The  PUV  gauges  provide  simultaneous  traces  of  submerged  pressure  p(f;x„,Zp)  and 
submerged  horizontal  velocity  components  iu,v)  or  u^(t-,x^,Zuy).  Both  instruments  are 
at  the  same  fixed  horizontal  position  (x,  y)  or  x„  but  at  different  (though  fixed  and 
adjacent)  vertical  elevations  Zp  and  z^y,  respectively,  as  sketched  in  Fig.  1.  For 
notation  convenience,  the  horizontal  tensor  notation  will  be  adopted;  the  position 
coordinates  are  (x^,z)  and  the  velocity  coordinates  (u^,w),  with  a=l  and  a  =  2 
corresponding  to  the  x  and  y  directions,  respectively,  of  a  cartesian  coordinate  system. 
The  z-axis  is  directed  vertically  upwards  from  the  horizontal  plane  of  the  local  mean 
water  level  (MWL),  in  opposition  to  the  gravity  vector  g. 

Short  extracts  from  two  measured  PUV  traces  are  shown  in  Fig.  2.  Trace  (a)  was 
measured  in  46  m  of  water  at  Platform  Edith,  about  9  km  offshore  from  Huntington 
Beach,  southern  California.  To  a  MWL  datum,  both  the  pressure  gauge  and  the 
directional  current  meter  are  at  —7.4  m.  The  sampling  rate  is  1  Hz.  Trace  (b)  was 
measured  in  17.6  m  of  water  off  the  Columbia  River  on  the  Pacific  coast  of  the  United 
States.  To  a  MWL  datum,  both  the  pressure  gauge  and  the  directional  current  meter  are 
at  — 16.8  m.  The  sampling  rate  is  4  Hz.  There  is  a  substantial  ambient  current  in  record 
(b). 

Although  combined  into  a  single  instrument,  the  precision  and  response  of  the 
pressure  and  directional  current  sensors  are  not  identical.  It  is  immediately  clear  from 
Lig.  2  that  there  is  significant  instrument  noise  in  the  UV  traces;  this  is  especially 


z  =  zp^ 


z  =  Zuv 


z  = 


-h 


Fig.  1.  Definition  sketch  and  PUV  gauge  location. 
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Fig.  2.  Measured  PUV  traces. 


notable  in  trace  (b)  where  the  sampling  interval  is  much  shorter.  The  P  traces  are 
consistently  smooth.  In  addition,  time  is  measured  independently  by  each  sensor,  and 
systematic  timing  discrepancies  are  possible.  This  is  not  evident  at  the  scale  of  these 
plots,  but  becomes  a  potential  concern  at  the  time  scale  of  the  local  analysis  that  is 
subsequently  outlined. 

Analysis  of  simultaneous  traces  such  as  Fig.  2  has  routinely  adopted  linear  wave 
theory  in  the  interpretation  of  the  local  kinematics  and  in  the  estimation  of  the  water 
surface  time  history.  This  prior  assumption  of  linearity  has  the  potential  for  significant 
misinterpretation  of  the  actual  local  kinematics.  Nonlinear  influences  may  be  filtered  and 
distorted  by  the  record  interpretation.  In  addition,  most  linear  analyses  assume  stationary 
conditions  over  a  measurement  span  of  sufficient  duration  to  obtain  adequate  frequency 
resolution  in  the  Fourier  transform. 

This  paper  will  introduce  and  demonstrate  a  nonlinear  theory  for  the  interpretation  of 
PUV  traces.  The  initial  discussion  will  review  and  illustrate  the  common  linear  analysis. 
This  leads  to  the  presentation  of  a  nonlinear  theory  for  the  interpretation  of  the  local 
irregular  kinematics  sampled  by  the  PUV  gauge.  The  theory  follows  the  spirit  of  the 
Sobey  (1992)  analysis  for  irregular  wave  kinematics  from  a  water  surface  record.  The 
nonlinear  PUV  theory  is  demonstrated  and  evaluated  for  theoretical  PUV  traces  from 
steady  wave  theory  and  for  measured  PUV  traces  from  the  Platform  Edith  and  Columbia 
River  sites  anticipated  in  Fig.  2. 
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2.  Global  linear  analysis 

Linear  analyses  of  direction-sensing  field  instruments  have  mostly  followed  the 
floating  buoy  analysis  of  Longuet-Higgins  et  al.  (1963)  (see  also  Horikawa,  1988).  This 
is  a  frequency  domain  analysis,  which  uses  the  entire  burst  sample  for  each  kinematic 
quantity  measured.  A  typical  burst  sample  has  a  duration  of  about  20  min,  such  that  it 
would  include  upwards  of  a  hundred  waves.  There  are  essentially  equivalent  implemen¬ 
tations  in  the  time  domain  (e.g.,  Lee  and  Wang,  1984). 

The  fundamental  basis  of  the  methodology  is  the  Gaussian  random  wave  model,  in 
which  the  irregular  water  surface  is  represented  as  the  superposition  of  very  many  linear 
waves  of  different  frequencies,  directions  and  amplitudes.  In  the  limit,  summing  over  all 
positive  frequencies  co  and  directions  6,  the  water  surface  rj  is  represented  as: 


1  ^17 

v(Aa,0  =  —  /  F(w,0)exp[-i(k„x„- wf)]dwd0. 

(1) 

where  Fia>,6)  is  the  complex  Fourier  transform  of  17  in  radian  frequency-direction 
space,  defined  by  the  inverse  Fourier  transform: 

^TT  ^00 

F(co,6-,x^,t)  =  1  I  r]{x^,t)&xp[i{k^x^—  cot)]dwd9, 

—  77*^  0 

(2) 

and  k^  are  the  cartesian  components  of  the  vector  wave  number  k.  k^  and 
through  the  linear  dispersion  relationship  for  waves  on  a  steady  current: 

CO  are  related 

(  w  -  k^U^Y  =  gk  tanh  kh. 

(3) 

in  which  is  the  local  depth-uniform  and  steady  current. 

It  follows  directly  from  linear  wave  theory  that  the  Fourier  transforms  of  the  dynamic 
pressure  and  the  horizontal  velocity  components  are: 

Fp{w,e-,x^,z,t)  =Kp{w,z)F{w,e-,x^,t), 

(4) 

=/:„(  w;z)yF(w,0;x„,f), 

(5) 

respectively,  where  k^/k  is  the  vector  (cos  6,  sin  9)  and: 

cosh  k(h  +  z)  cosh  k(h  +  z) 

(6) 

The  Kp  and  transfer  functions  are  not  dependent  on  the  observations.  In  practice, 
the  Doppler  adjustment  for  local  current  is  mostly  omitted,  a  common  practice  being  to 
remove  the  mean  level  (i.e.,  the  current)  from  the  velocity  traces  prior  to  analysis. 

Longuet-Higgins  et  al.  (1963)  approximated  the  local  directional  variance  spectrum 
Eiw,0),  strictly  E^^(w,9),  as  the  first  five  terms  in  a  Fourier  series  in  6: 

2  1 

E(w,6)  =  Oq  H —  (  cos  0  -I-  sin  0)  H —  (02  cos 20  -I-  f>2  sin20)  -I-  . . .  (7) 
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in  which  (Grosskopf  et  ah,  1983;  Horikawa,  1988)  the  coefficients  are  frequency-depen- 
dent: 


2o{co)  = 


2l(  W)  = 


hi(w)  = 


Eppi  «) 


2ttKI{w,Zp) 


Ep3  ") 


tt/:  (  w;zp)/:„(w;zuv) 


Epul ") 


77/:  (w;zp)/:„(w;zuv) 


2(")  = 


b^ia)  = 


Eu^ul^  tu) 

a);zuv) 


(8) 

(9) 

(10) 


£pp((u),  £„  „(&))  and  £„  „  (&))  are  the  auto-variance  spectra  estimated  from  the  p,  Mj 
and  M2  burst  samples.  Ep^j^co),  (w)and  „  (cu)  are  the  cross-variance  spectra  from 
the  same  burst  samples. 

This  algorithm  has  the  advantages  of  convenience,  relative  simplicity  and  many  years 
of  routine  use.  The  fundamental  reliance  on  linear  theory  and  on  spectral  analysis  are 
both  potentially  significant  constraints.  Linear  theory  compromises  the  identification  and 
extraction  of  nonlinear  influences  in  the  data.  The  practice  of  spectral  analysis  insists 
that  the  local  kinematics  at  time  t  be  dependent  with  almost  equal  weight  on  all  data 
observations  at  t  =  0,  Af,2At,3At, . . .  ,{N  —  l)Af  in  the  burst  sample.  For  PUV  data,  Af 
is  typically  1  s  and  the  NAt  duration  of  the  data  is  typically  20  min,  for  which  N 
exceeds  1000  observations.  This  is  a  linear  and  very  global  (burst-sample  duration) 
interpretation.  To  preserve  the  physical  integrity  of  the  data,  a  very  local  interpretation 
would  be  preferable.  Locally  linear  and  then  locally  nonlinear  theories  are  presented  in 
the  following  sections. 


3.  A  local  linear  analysis 


A  linear  but  locally  focused  interpretation  is  a  useful  intermediate  step,  in  the  spirit  of 
the  Nielsen  (1989)  locally  linear  approximation  to  irregular  wave  kinematics  from  a 
water  surface  trace.  The  basis  is  linear  theory  estimates  of  dynamic  pressure  and 
horizontal  velocity  components: 

Pa{x^,z. 

cosh  k{h  +  z) 

coshtt 

(11) 

u^{x„,z 

k„  cosh  k(h  +  z) 

.,)-[/,+  t.c/j  ,(*,.<). 

(12) 

in  which: 

V(x^,t) 

=  a  cos(  k^x^  —  wt). 

(13) 

is  the  linear  water  surface,  and  k  and  w  are  related  through  the  linear  dispersion 
relationship  (Eq.  (3)). 
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A  prediction  of  the  local  kinematics  is  sought  throughout  the  water  column  in  the 
immediate  neighborhood  of  the  PUV  gauge  at  horizontal  position  x^.  The  pressure 
sensor  is  at  known  elevation  Zp  and  the  directional  current  meter  at  known  elevation 
Zjjv-  Measured  dynamic  pressure  traces  =  Pd(trx^,Zp)  are  available  at  discrete 
times  f;,  as  are  measured  velocity  component  traces  m*®  =  uj^tf,x^,z^y)  at  the  same 
discrete  times.  The  local  water  depth  h  and  depth-uniform  current  are  also  known. 
The  current  can  be  estimated  from  the  UV  traces  as  a  time-average  over  sufficient  time 
to  average  over  the  local  waves  but  not  over  local  astronomical  and  storm  tide  currents. 

At  each  time,  there  are  four  unknowns  ico,k^  and  rj)  and  four  equations,  dispersion 
(Eq.  (3))  plus  the  three  observational  equations  (Eqs.  (ll)  and  (12)).  These  are  the 
simultaneous  implicit  algebraic  equations: 


/i(w,k„,i?)  =(w-k„t4)  -gktanhkh  =  Q 
cosh  k(h  +  Zp) 

/,( -  rt- -  0 

k.  cosh  k(/i -I- Zttv) 


(14) 


f,{w,k^,r])  =  -{oi-k^U^) 


cosh  k{h  +  Zuv) 
sinh  kh 


V  -  wf  *  =  0, 


In  principle,  a  unique  local  solution  at  each  time  f,  can  be  computed.  An  analytical 
solution  does  not  seem  feasible.  A  direct  numerical  algorithm  is  the  Newton-Raphson 
method.  Writing  Eq.  (14)  as  fi(xj)  =  0  where  i,j  =  1, 2,3,4  and  vector  as  the 
solution  estimate  at  iteration  n,  the  correction  is  suggested  by  the  local  Taylor 

series  expansion: 


df^ 

dxj 


(«) 

...  =0. 


(15) 


The  Jacobian  Ofj/Oxj  is  evaluated  analytically  and  Eq.  (15)  is  solved  by  matrix 
inversion. 

The  predictive  potential  of  this  algorithm  was  evaluated  for  a  theoretical  steady  wave 
train  of  period  10  s  and  height  10  m,  directed  at  -I-  tt/IO  to  the  x-axis  in  20  m  of  water 
and  an  opposing  current  of  —  1  m/s.  The  PUV  gauge  is  located  at  Zp  =  Zyy  =  —10  m. 
The  initial  theoretical  PUV  trace  was  computed  from  Airy  theory  at  a  sampling  time 
interval  of  0.5  s.  With  this  trace  as  the  measured  PUV  record,  the  local  linear  algorithm 
predicted  the  local  wave  frequency,  the  local  wave  number  components  and  the  local 
water  surface  elevation.  Except  in  the  immediate  neighborhood  of  zero-crossings,  where 
profile  curvature  is  minimal  and  Eq.  (14)  is  ill-conditioned,  there  was  consistent 
agreement  with  Airy  theory  to  four  significant  figures. 

A  second  theoretical  PUV  trace  (record  ‘Twenty’;  see  Table  1)  was  computed  from 
near-exact  (global)  Eourier  wave  theory  (Sobey,  1989),  also  at  a  sampling  time  interval 
of  0.5  s.  With  this  trace  as  the  measured  PUV  record,  the  local  linear  algorithm  was 
significantly  less  successful. 
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Table  1 

Sample  theoretical  waves 


Wave 

H 

h 

T 

U, 

Zp 

Zuv 

Five 

3  m 

5  m 

10  s 

-1.902 

-0.618 

—  5.0  m 

—  5.0  m 

Twenty 

10  m 

20  m 

10  s 

-0.951 

-0.309 

- 10.0  m 

- 10.0  m 

Hundred 

20  m 

100  m 

10  s 

-0.478 

-0.155 

-20.0  m 

-20.0  m 

The  predicted  water  surface  elevation  and  wave  number  components  from  the  local 
linear  algorithm  are  shown  in  Fig.  3  as  the  markers,  together  with  the  near-exact 
predictions  from  global  Fourier  wave  theory  as  the  solid  lines.  Solutions  to  Eq.  (14)  at 
times  in  the  immediate  neighborhood  of  the  zero-crossings  were  difficult  to  obtain,  and 
required  initial  solution  estimates  that  were  almost  exact.  For  this  ‘wave  transitional’ 
situation,  a  solution  was  obtained  throughout.  It  is  smooth  and  visually  convincing, 
except  when  compared  with  the  near-exact  predictions.  There  is  order  of  magnitude 
agreement  only.  The  crest  elevation  in  particular  is  poorly  predicted. 

Some  variations  on  the  locally  linear  formulation  may  achieve  a  more  acceptable  and 
perhaps  also  more  robust  local  solution,  but  the  indications  are  not  encouraging. 
Eliminating  17  between  the  /j  and  both  the  /j  and  /4  equations  gives: 


t(s) 


Fig.  3.  Local  linear  theory  (markers)  and  global  Fourier  theory  (solid  lines)  predictions  for  Record  ‘Twenty’. 
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Squaring  now  and  using  in  the  /j  equation  establishes  a  single  implicit  algebraic 
equation  in  the  local  wave  number  k: 


pg  cosh  k{h  +  Zp) 
cosh  k{h  +  Zuv ) 


gk 


tanh  kh 


=  0.  (17) 


Some  special  cases  of  this  equation  are  very  revealing.  When  Zp  =  z^y,  Eq.  (17) 
becomes: 


pg 


kh 


tanh  kh 


0. 


(18) 


The  function  kh /tanh  kh  is  1  at  kh  =  0  (‘wave  shallow’  water),  and  approaches  kh 
asymptotically  from  above  for  large  kh  (‘wave  deep’  water).  There  are  solutions  only 
for: 


1 

gh 


pg 


>  1. 


(19) 


Solutions  do  not  exist  under  quite  a  wide  range  of  conditions.  This  is  a  potentially  fatal 
impediment. 

It  may  be  possible  to  reformulate  the  locally  linear  problem  to  avoid  these  difficul¬ 
ties.  Those  time  steps  where  there  is  no  solution  might  be  accommodated  by  the  solution 
of  Eq.  (14)  in  the  least-squares  sense  rather  than  an  exact  solution.  It  may  also  be 
advisable  to  use  neighboring  observations  in  the  local  solutions.  Erom  the  experience  of 
the  Sobey  (1992)  locally  nonlinear  solution  from  water  surface  traces,  both  of  these 
measures  might  help  in  dealing  with  observational  error  bands.  All  of  these,  however, 
increase  the  complexity  of  the  locally  linear  formulation  and  negate  the  sole  advantage 
of  a  linear  formulation,  simplicity.  A  locally  nonlinear  analysis  is  viable  and  is 
introduced  in  Section  4. 


4.  A  local  nonlinear  analysis 


Nonlinear  and  irregular  waves  follow  a  field  equation,  the  Laplace  equation: 
8^4, 

- 1 - g-  =  0, 

dx^dx^  dz 

together  with  kinematic  and  dynamic  free  surface  boundary  conditions: 
drj  drj 

fK  =  W-  Z  =  V, 

ot  dx^ 

/d  =  —  +  -(MpMp  +  >v2)  +  gr]  =  Oat  Z=  V, 

and  a  kinematic  bottom  boundary  condition  on  a  sloping  bed: 
dh 


w  +  u 


“dx^ 


'  Oat  z  =  —h. 


(20) 

(21) 

(22) 


(23) 
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These  equations  assume  only  that  the  flow  is  irrotational  and  incompressible.  The 
velocity  potential  function  is  defined  such  that  the  velocity  components  are  (u^,w)  = 
(9(/)/9x„,9(/)/8z). 

If  it  is  further  assumed  that  the  bottom  is  locally  horizontal,  the  bottom  boundary 
condition  becomes  just: 

w  =  Oat  z  =  —h.  (24) 

Nonlinear  wave  solutions  of  the  form: 

cosh  jk{h  +  z) 

=  U^x^  +  l,Aj - - &iT^j{kaX^-wt),  (25) 

exactly  satisfy  the  field  (Eq.  (20))  and  the  reduced  bottom  boundary  condition  (Eq. 
(24)).  The  current  must  be  steady  and  depth-uniform,  excluding  any  consideration  of 
velocity  shear. 

Erom  the  known  velocity  potential  function,  the  velocity  components  are: 

cosh  jk{h  +  z) 

=  U^  +  Uk^A. - - cos  j{k^x^-  wt),  (26) 

sinh  jk(  h  +  z) 

w{x„,z,t)  =IijkAj - - sin  j{k^x^-  wt).  (27) 

Using  also  the  irrotational  Bernoulli  equation,  the  dynamic  pressure  is: 

_(?</)  1 

—  +  (28) 

in  which: 


dej)  sinh  jk(h  +  z) 

^  =  -  Di - TTUTiT - X  k^x^-wt). 


cosh  jkh 


and  the  Bernoulli  constant  is  (Sobey,  1992): 


jkAj 


cosh  jkh  j 


(29) 


(30) 


A  prediction  of  the  local  kinematics  is  sought  throughout  the  water  column  in  the 
immediate  neighborhood  of  the  PUV  gauge  at  horizontal  position  x^.  The  pressure 
sensor  is  at  known  elevation  Zp  and  the  directional  current  meter  at  known  elevation 
Zuv  Measured  dynamic  pressure  traces  =  Pi(ti-,x^,Zp)  are  available  at  discrete 
times  f;,  and  measured  velocity  component  traces  m*®  =  M^(t;;x„,Zuv)  ^^e  available  at 
the  same  discrete  times.  The  local  water  depth  h  and  depth-uniform  current  U^,  are  also 
known.  The  current  can  be  estimated  from  the  UV  traces  as  a  time-average  over 
sufficient  time  to  average  over  the  local  waves  but  not  over  local  astronomical  and  storm 
tide  currents. 

The  unknowns  in  this  local  solution  are  radian  frequency  w,  the  wave  number 
components  k^,  the  spatial  phase  k^x^,  the  Eourier  coefficients  A^  and  the  local  water 
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surface  elevations  fln~  neighborhood  of  time  f,.  With  J  Fourier 

coefficients,  N  local  water  surface  elevations,  and  treating  the  spatial  phase  as  a  single 
unknown,  there  are  A  +  J  +  N  unknowns. 

Active  equations  are  provided  by  the  PUV  observations  and  by  the  kinematic  and 
dynamic  free  surface  boundary  conditions.  Note  that  the  primitive  form  of  the  kinematic 
free  surface  boundary  condition,  Eq.  (21),  requires  estimates  of  drj/dt  and  drj/dx^.  The 
temporal  gradient  'drj/dt  might  be  estimated  from  interpolation  among  the  local 
which  is  part  of  the  solution,  but  there  is  no  spatial  information  from  the  PUV  traces. 
Sobey  (1992)  estimated  these  spatial  gradients  by  imposing  a  locally  steady  approxima¬ 
tion  on  the  waves. 

This  additional  approximation  can  be  avoided.  The  problem  formulation  is  not 
compromised  by  redefinition  (Longuet-Higgins,  1962)  of  the  kinematic  free  surface 
boundary  condition  as: 


/k  /k 


1  Df^ 


8  Dt 


1  D 
=  w  H - 


I  d(j}\  Mp  Z)Mp  w  Dw 

I  I  _  r\^  _  r\^ 


(31) 


This  modified  form  of  the  kinematic  free  surface  boundary  condition  excludes  both 
temporal  and  spatial  gradients  of  p. 

As  an  initial  evaluation  of  this  form  of  the  kinematic  free  surface  boundary  condition, 
it  was  substituted  in  the  Sobey  (1992)  code  for  local  irregular  wave  kinematics  from  a 
single  water  surface  trace  (local  Fourier  irregular;  LFI-E).  For  theoretical  traces  where 
the  locally  steady  approximation  was  exact,  the  results  were  visually  identical.  For 
measured  traces,  the  results  were  similar.  Without  the  need  to  estimate  drj/dt  from 
measured  water  surface  elevations,  the  local  numerical  solutions  appeared  to  be  much 
more  robust.  The  Eq.  (31)  form  of  the  kinematic  free  surface  boundary  condition  was 
adopted  in  the  present  PUV  analysis. 

At  any  time  t  where  there  is  an  assigned  water  surface  node,  active  theoretical 
equations  are  provided  by  both  free  surface  boundary  conditions.  In  addition,  there  are 
PUV  observational  equations  at  each  discrete  measurement  time  f,.  Each  of  these 
equations  are  implicit  and  algebraic.  Available  equations  are  the  modified  KFSBC  at  rj, 
the  DFSBC  at  17,  a  pressure  equation  at  Zp  and  two  horizontal  velocity  equations  at 
Zuv  • 

I  D  I  d(l)\  u„  Du„  w  Dw 

/K(w,k„,k„x„,A.,T7)  =w-f  -  —  —  -f - —  -f  -  —  =0at  z=  7? 

g  Dt  \  dt  j  g  Dt  g  Dt 

1  / 

fp,{w,k^,k^x^,Aj,r])  =  —  +  -{ul  +  w^)+g7]  =  0atz=V 
fp(co,k^,k^  x^,Aj,Tq)  =p„-p*‘‘  =  0at  z  =  zp 
/J ^^k^,k^x^,A.,r])  =  u^-  =  0 at  z  =  Zuv 


(32) 
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Additional  closure  equations  can  be  provided  by  observational  and  free  surface  bound¬ 
ary  equations  at  neighboring  times  within  a  short  local  window  of  duration  Tq  that  is 
centered  on  f,.  At  least  4  +  J  +  N  independent  equations  must  be  provided  to  solve  for 
4  +  J  +  N  unknowns. 

Apart  from  this  strict  mathematical  closure  requirement,  there  is  an  additional 
constraint.  The  K  and  D  free  surface  boundary  equations  are  exact,  but  the  P,  U  and  V 
observational  equations  have  measurement  error  bands.  Also,  the  P  and  U,  V  sensors  are 
unlikely  to  have  the  same  accuracy.  Present  technology  has  much  larger  error  bands  on 
the  U,  V  traces  than  on  the  P  traces;  this  is  seen  quite  clearly  in  Fig.  2b,  where  the 
sampling  rate  is  4  Hz.  These  realities  of  the  problem  formulation  are  accommodated  by 
a  significant  overspecification  of  the  problem,  especially  with  the  UV  observational 
equations,  and  the  adoption  of  a  least  squares  rather  than  an  exact  solution.  Flexibility  in 
the  time  location  of  the  PUV  observational  equations  was  also  introduced  through  cubic 
spline  interpolation  among  the  measured  points. 

The  local  LFI-PUV  theory  has  three  free  parameters,  the  truncation  order  J  and  the 
number  of  17  points  N  in  each  local  window  together  with  the  width  Tg/7i,  of  the  local 
windows.  The  truncation  order  has  much  the  same  authority  as  order  in  an  analytical 
wave  theory  (Stokes,  conoidal)  or  truncation  order  in  Fourier  wave  theory.  A  window 
width  of  Tg/T^  =  1  would  be  a  global  solution.  A  particular  solution  will  be  designated 
LFI-PUV[7,A,To/7;].  For  example,  LFI-PUV[3,3,0.1]  has  7=3,  N=3  and  Tq/T;  = 
0.1. 


5.  Numerical  implementation 

As  in  Sobey  (1992),  the  analysis  segment  routinely  adopted  was  a  double  wave 
sequence  centered  about  a  crest  suggested  by  the  pressure  trace.  Solutions  are  sought  in 
narrow  local  windows,  with  a  target  width  of  order  Tq  =  O.IT,,  where  T,  is  the  zero-up 
crossing  period  of  the  record  segment. 

In  all  applications  of  the  present  theory,  it  has  proven  convenient  to  assign  an  odd 
number  of  local  water  surface  elevation  i7(f„);  with  N=1  or  3  (for  a  single  width 
window)  or  5  (for  a  double  width  window),  such  that  there  is  always  a  computed  water 
surface  elevation  centrally  located  in  the  window.  With  respect  to  the  time  at  the  center 
of  a  window  of  width  Tq,  rjitJ  points  are  located  at  J/tq  =  0;  +0.5;  +1.0.  K  and  D 
equations  are  applied  at  each  of  these  rjitj  points,  depending  on  the  value  of  N.  At 
N  =  5,  the  local  window  width  is  doubled.  There  are  2N=2  or  6  (for  a  single-width 
window)  or  10  (for  a  double-width  window)  of  these  water  surface  equations. 

P,  U  and  V  observational  equations  are  located  at  the  center  of  the  window. 
Additional  PUV  equations  within  the  window  are  located  such  that  UV  observational 
equations  are  not  assigned  a  weighting  that  was  inconsistent  with  their  routinely  larger 
error  bands.  With  respect  to  the  time  at  the  center  of  a  window  of  width  Tq,  P  equations 
are  located  at  f/r,,  =  +0.25,  +0.5  in  a  single  width  window  and  also  at  +0.75,  +  1  in 
a  double  width  window.  The  associated  U  and  V  equations  are  located  at  J/tq  =  +0.125, 
+  0.25,  +0.375,  +0.5  in  a  single-width  window  and  also  +0.625,  +0.75,  +0.875, 
+ 1  in  a  double-width  window.  The  explicit  higher  density  of  UV  equations,  together 


28 


RJ.  Sobey,  S.A.  Hughes  /  Coastal  Engineering  36  (1999)  17-36 


with  the  least  squares  algorithm,  was  adopted  as  partial  accommodation  of  the  larger 
error  bands  of  the  UV  observational  equations.  The  least  squares  objective  function  was 
normalized  as: 


0{(o,k„,k^x„-,Afrin)  =  ^S(/k+/d)  +  +  ^^^(/u+/v). 

(33) 

such  that  K,  D,  P  and  UV  equation  groups  have  equal  weight.  Np  is  1  +  4L  and  N^y  is 
1  +  8L,  where  L  is  1  (single-width  window)  or  2  (double-width  window). 

There  are  3  -I-  20  L  =  23  (for  a  single  width  window)  or  43  (for  a  double  width 
window)  observational  equations.  There  are  a  total  of  2V-I-3-I-20L  physical  and 
observational  equations  in  4  +  J  +  N  unknowns,  where  2N  +  3  +  2QL  must  equal  or 
exceed  4  +  J  +  N.  A  typical  application  has  7=3  and  N  =  3,  so  that  there  are  10 
unknowns:  Aj, A3, pj, 172,173.  For  a  single-width  window,  there  are  29 

equations,  and  49  for  a  double-width  window.  Experience  suggested  that  this  level  of 
overspecification  was  generally  necessary  to  accommodate  the  observational  error 
bands.  For  theoretical  wave  traces  with  no  observational  error  bands  (see  Section  6), 
only  a  closed  set  of  KD  and  PUV  equations  was  necessary  for  credible  solutions; 
overspecification  did  not  compromise  this  success. 

The  problem  formulation  requires  the  least-squares  solution  of  a  system  of  simultane¬ 
ous,  nonlinear,  implicit  algebraic  equations.  Numerical  solutions  routinely  used  the 
Dennis  et  al.  (1981b)  NL2SOL  code;  an  updated  code  is  available  from  the  internet 
Netlib  repository  of  scientific  subroutines  in  fortran.  The  algorithm  (Dennis  et  al., 
1981a)  is  a  variation  on  Newton’s  method  in  which  part  of  the  Hessian  matrix  is 
computed  exactly  and  part  is  approximated  by  a  secant  updating  method.  To  promote 
convergence  from  poor  starting  guesses,  a  model/trust-region  technique  is  used  along 
with  an  adaptive  choice  of  the  model  Hessian.  In  operation,  the  algorithm  sometimes 
reduces  to  the  classical  Gauss-Newton  or  Levenberg-Marquardt  methods. 

The  number  of  unknowns,  typically  10,  is  relatively  large.  In  any  such  high  order 
nonlinear  optimization  problem,  the  crucial  elements  of  successful  solutions  are  invari¬ 
ably  accurate  estimates  of  the  Jacobian  and  good  initial  solution  estimates.  Accurate 
estimates  of  the  Jacobian  require  analytical  estimates  of  the  partial  derivatives  of  all  / 
equations  with  respect  to  all  the  unknown  parameters.  These  are  the  quantities 
df^/d(o,df^/dki,df^/dk2,  ■  ■  ■  through  ...,  df^^/dr]^^.  The  analytical  derivatives  were 
evaluated  and  confirmed  against  finite  difference  approximations. 

Establishing  good  initial  solution  estimates  always  requires  extensive  experience  with 
the  particular  physical  problem  and  system  of  equations.  The  methodology  finally 
adopted  was  a  two-stage  algorithm.  The  initial  step  identified  three  overlapping  record 
sub-segments,  the  half-wave  extract  centered  on  the  leading  trough,  the  half-wave 
centered  on  the  central  crest  and  the  half-wave  centered  on  the  following  trough.  Each 
half-wave  is  equivalent  to  an  extra-wide  window  of  width  about  O.SIi,.  In  each 
half-wave  of,  say,  yk"  discrete  PUV  observations,  17  points  were  located  at  the  same 
times  such  that  there  are  4  -I-  7  -I-  ,7^  unknowns.  K,  D,  P,  U  and  V  equations  are  assigned 
at  the  same  times,  giving  5yF  equations.  For  a  PUV  record  sampled  at  3  Hz,  there 
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would  be  about  observational  points  in  a  half-window.  For  7=3  and  3=1,  there 
would  be  14  unknowns  and  35  equations.  For  7=3  and  3=4,  there  would  be  35 
unknowns  and  140  equations.  With  initial  estimates  provided  by  a  global  Airy  theory 
approximation,  smooth  and  robust  solutions  are  established  in  each  extra-wide  half¬ 
window.  These  in  turn  provide  routinely  appropriate  initial  estimates  for  the  local 
solutions  in  very  much  narrower  windows. 

The  algorithm  had  most  difficulty  in  finding  credible  local  solutions  around  the 
zero-crossings  of  the  P  trace.  At  such  regions,  profile  curvature  is  very  small  and  a 
narrow  local  window  has  poor  resolution.  Often,  this  is  not  a  region  of  particular 
concern.  Two  strategies  are  possible  in  dealing  with  this  problem.  The  first  would  be  to 
ignore  such  points  and  rely  on  interpolating  for  the  kinematics  from  adjacent  local 
solutions  that  bracket  the  zero-crossings.  The  second  would  be  to  extend  the  local 
window  to  improve  resolution  of  the  local  profile  curvature.  This  second  approach  has 
been  adopted,  and  was  the  essential  rationale  for  the  double-width  window  introduced 
above. 


6.  Theoretical  wave  traces 

An  initial  evaluation  of  the  LFI-PUV  theory  and  coding  is  provided  by  theoretical 
PUV  traces  from  uniform,  long-crested  wave  trains.  Three  monochromatic  waves. 
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‘Five’,  ‘Twenty’  and  ‘Hundred’,  have  been  defined,  respectively,  in  wave  shallow,  wave 
transitional  and  wave  deep  water.  The  details  are  listed  in  Table  1.  The  wave  heights  in 
all  three  cases  are  moderately  extreme  for  the  depth.  Note  in  particular  the  location  of 
both  Zp  and  Zuv  higher  in  the  water  column  for  waves  ‘Twenty’  and  ‘Hundred’.  The 
significant  attenuation  in  wave  kinematics  with  depth  is  frequently  mitigated  by  such  a 
practice.  Global  Fourier  wave  theory  (Sobey,  1989),  at  truncation  order  18  with  25  water 
surface  nodes  from  crest  to  trough,  provides  near-exact  kinematics  for  these  wave  trains. 

The  PUV  observational  traces  at  the  listed  elevations  were  predicted  at  a  time  interval 
of  0.5  s.  LFI-PUV[3, 3,0.1]  predictions  for  wave  ‘Twenty’  are  shown  in  Fig.  4.  Part  (a) 
shows  the  ‘measured’  P(|(f;Zp)  trace  together  with  the  ‘measured’  Mj(f;zuy)  and 
U2(t;zi]y)  traces.  Part  (b)  shows  the  rj  traces  predicted  by  LFI-PUV  (markers)  and  by 
the  near-exact  Fourier  wave  theory  (continuous  line).  Compare  this  result  to  the 
relatively  poor  prediction  of  17  for  wave  ‘Twenty’  using  the  locally  linear  method  (Fig. 
3).  Part  (c)  shows  all  three  velocity  components  at  elevation  z=  —H/2  just  below  the 
trough,  and  part  (d)  all  three  acceleration  components  also  at  the  same  near  trough 
elevation.  Again,  the  markers  are  LFI-PUV  predictions  and  the  continuous  lines  are 
from  the  near-exact  Fourier  wave  theory.  Agreement  throughout  is  almost  perfect. 

The  excellent  agreement  is  achieved  at  a  relatively  low  order  (7=3),  and  this  is  a 
consistent  strength  of  the  LFI  methodology.  This  success  at  such  low  order  is  achieved 
by  seeking  separate  solutions  in  each  local  window.  The  solution  parameters, 
will  vary  from  window  to  window.  The  variation  in  the  w,  k^  and 
Aj  parameters,  that  would  also  be  defined  in  a  global  solution,  is  not  expected  to  be  very 
large.  The  local  variation  in  these  solution  parameters  for  wave  ‘Twenty’  is  shown  in 
Fig.  5.  Note  the  double-width  windows  at  the  profile  zero-crossings  near  t=  +2  s. 
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Fig.  5.  LFI-PUV[3,3,0.1]  local  solution  parameters  for  wave  ‘Twenty’. 
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Fig.  6.  LFI-PUV[3,3,0.1]  predictions  for  wave  ‘Five’. 


Similarly  excellent  agreement  is  achieved  for  waves  ‘Five’  and  ‘Hundred’,  respec- 
lively,  in  wave  shallow  and  wave  deep  water.  Comparisons  of  LFI-PUV[3,3,0.1] 
solutions  and  near-exact  Fourier  kinematics  are  shown  in  Figs.  6  and  7.  Collectively, 
waves  ‘Five’,  ‘Twenty’  and  ‘Hundred’  demonstrate  the  suitability  of  the  LFI-PUV 
algorithm  under  expected  operational  conditions  from  wave  shallow  through  wave 
transitional  to  wave  deep  water. 


7.  Measured  wave  traces 

The  theoretical  traces  of  Section  6  provide  confidence  in  the  problem  formulation  and 
in  the  numerical  algorithm.  But  the  LFI-PUV  methodology  is  intended  for  the  nonlinear 
interpretation  of  field  PUV  traces.  The  inherent  measurement  error  bands  of  field  PUV 
traces  will  provide  a  much  more  challenging  test. 

The  two  field  records  partially  illustrated  in  Fig.  2  have  been  used  for  evaluation. 
Apart  from  sampling  very  different  wave  climates,  these  records  differ  very  significantly 
in  the  sampling  rates,  1  Hz  for  the  Platform  Edith  data  and  4  Hz  for  the  Columbia  River 
data. 

Figs.  8  and  9  show  sample  LFI-PUV[3,3,0.2]  predictions  of  local  kinematics  from 
these  records.  These  predictions  are  credible  and  they  match  the  measured  PUV  traces 
by  definition,  but  there  is  no  independent  field  data  for  confirmation.  In  the  Sobey 
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Fig.  8.  LFI-PUV[3,3,0.2]  predictions  from  Platform  Edith  data. 
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^(b)  Velocity  components  at  ^  =  —2  m 


(1992)  LFI-E  theory  for  irregular  wave  kinematics  from  single-point  water  surface 
traces,  a  window  width  of  Tq/C,  =  0.1  was  routinely  satisfactory.  For  PUV  traces,  the 
observational  equations  are  somewhat  less  authoritative.  An  obvious  reason  is  the 
relatively  large  error  bands  on  the  UV  measurements.  Another  possible  reason  is  natural 
signal  attenuation  with  depth  of  submergence;  the  higher  frequencies  are  most  rapidly 
attenuated  and  their  influence  may  get  lost  in  the  error  bands  (Barker,  1998).  Neverthe¬ 
less,  experience  suggested  the  routine  adoption  of  a  slightly  wider  window  width  of 
Tq/T^  =  0.2.  Barker  and  Sobey  (1996)  also  recommended  this  same  window  width  for 
irregular  wave  kinematics  from  pressure  records  alone. 

It  is  instructive,  nevertheless,  to  compare  the  locally  nonlinear  prediction  of  sea 
surface  elevation  with  the  locally  linear  prediction  (shown  as  the  markers  on  Figs.  8b 
and  9b).  The  density  of  the  linear  solutions  match  the  data  resolution.  In  Fig.  8b,  Af  is 
1.0  s.  Where  there  is  no  solution  (see  Fq.  (17)  and  following),  the  marker  is  located  at 
the  MWF,  as  at  times  —4,  —3  —2,  2,  6,  7  and  10  s.  Solutions  at  0,  1  and  9  s  are 
significantly  off  scale,  and  only  those  at  —  1,  3,  4,  5  and  8  s  have  any  credibility.  In  Fig. 
9b,  At  is  0.25  s.  There  are  many  ‘no  solution’  markers  and  the  few  locations  where 
linear  solutions  were  obtained  are  mostly  off-scale  or  have  little  credibility.  The 
theoretical  reservation  on  the  linear  method,  potentially  significant  underprediction  at  the 
crest  (see  Fig.  3),  is  clearly  overwhelmed  by  the  frequency  of  either  ‘no  solution’  results 
or  solutions  that  have  no  visual  credibility.  The  linear  method  must  be  rejected  on  both 
theoretical  and  practical  grounds. 
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The  LFI-PUV  solutions  for  the  Columbia  River  data  have  a  somewhat  tenuous 
nature,  evidenced  by  the  sharp  changes  in  predicted  kinematics  in  the  neighborhood  of 
times  —8,  +5  and  +8  s.  The  UV  error  bands  for  this  record  (Fig.  2b)  are  quite 
significant.  The  specific  PUV  gauge  was  designed  with  the  expectation  that  analysis 
would  follow  the  Longuet-Higgins  et  al.  (1963)  linear  and  global  frequency-domain 
analysis  outlined  in  Eqs.  (I)-(IO).  In  this  global  linear  analysis,  measurement  error  is 
easily  accommodated  by  spectral  smoothing  of  the  £^^(^(0,6)  through  E^^^^co,0) 
auto-variance  and  cross-variance  spectra  and/or  the  E^^io),0)  estimate  from  Eq.  (7). 

Measurement  error  bands  are  not  so  easily  dismissed  in  the  local  LEI-PUV  theory. 
LEI-PUV  seeks  an  enhanced  representation  of  the  local  kinematics,  through  the  direct 
use  of  local  PUV  measurements  in  the  problem  formulation.  The  local  PUV  measure¬ 
ments  provide  the  local  reality,  so  that  the  predictive  potential  can  be  significantly 
compromised  by  local  measurement  error.  The  Columbia  River  data  set  provides  a 
useful  reality  check.  While  the  value  of  the  Pig.  9  predictions  are  arguably  tenuous,  their 
particular  value  may  lie  in  the  identification  of  UV  measurement  error  bands  as  a 
relevant  design  issue  for  PUV  instrumentation. 

An  associated  instrument  design  issue  is  the  attenuation  of  the  fluctuating  kinematics 
(P  and  UV)  with  increasing  depth  of  submergence  and  with  increasing  frequency.  The 
attenuation  increases  with  both  frequency  and  depth  of  submergence.  Transducer 
accuracy  is  absolute,  so  that  relative  accuracy  is  significantly  compromised  for  deeply 
submerged  instruments  at  higher  frequencies  in  deeper  water.  There  is  only  partial 
mitigation  for  this  in  the  superior  accuracy  of  P  measurements,  with  respect  to  UV 
measurements.  The  design  compromise  has  been  to  locate  PUV  gauges  reasonably  high 
in  the  water  column  to  avoid  this  problem  as  much  as  possible.  The  Columbia  River 
gauge  is  deployed  at  —  16.8  m  where  the  water  depth  is  17.6  m,  but  the  Platform  Edith 
gauge  is  deployed  at  only  —7.4  m  where  the  water  depth  is  46  m.  The  deeper 
submergence  of  the  Columbia  River  gauge  contributes  to  the  tenuous  nature  of  the 
analysis  of  this  data. 

A  careful  review  of  the  Pig.  9  predictions  reveals  that  the  maximum  vertical  velocity 
prediction,  1.59  m/s  at  -1-0.5  s,  exceeds  the  maximum  horizontal  velocity  prediction, 
1.39  m/s  at  —0.5  s.  In  a  global  steady  wave  solution,  such  as  the  ‘global  Pourier’ 
solutions  in  Pigs.  4,  6  and  7,  the  maximum  horizontal  and  vertical  velocities  will  again 
not  coincide  in  time,  but  the  maximum  vertical  velocity  magnitude  will  always  be 
smaller  than  the  maximum  horizontal  velocity  magnitude.  This  underlines  the  local 
character  of  the  PPI-PUV  methodology.  The  PPI-PUV  solutions  in  each  local  window 
are  independent,  and  directly  reflect  the  local  measured  kinematics.  There  is  no  global 
influence. 

A  final  phase  in  evaluation  of  the  PPI-PUV  methodology  would  be  detailed  labora¬ 
tory  measurements  of  kinematics  under  regular  and  irregular  wave  conditions.  These 
should  include  a  PUV  instrument  plus  the  local  water  surface  and  (m„,w)  velocity 
component  and  dynamic  pressure  measurements  at  one  or  more  different  near-surface 
elevations.  Suitable  laboratory  data  do  not  presently  exist.  Establishing  a  suitable  data 
set  would  be  a  very  major  task.  It  would  presumably  involve  sophisticated  three-dimen¬ 
sional  laser  Doppler  instrumentation  for  the  velocity  measurements.  Pressure  is  very 
rarely  measured  in  a  wave  laboratory,  except  at  structures  where  the  data  would  be 
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unsuitable.  Pressure  measurement  at  near-surface  elevations  under  waves  would  require 
special  attention. 


8.  Conclusions 

Traditional  directional  analysis  of  PUV  data  results  in  linearized  estimates  of 
directional  wave  and  velocity  spectra  which  typically  do  not  include  steady  currents. 
Attempts  to  reconstitute  velocity  time  series  at  any  other  location  in  the  water  column 
result  in  linear  estimates  with  decreasing  veracity  approaching  the  free  surface.  The 
single  advantage  of  linear  directional  spectral  analysis  is  the  compactness  of  the  output 
which  corresponds  to  loss  of  information.  However,  this  compactness  could  become  a 
liability  in  some  applications. 

A  new  method  has  been  presented  for  interpreting  traces  recorded  by  conventional 
PUV  instruments  under  irregular  multidirectional  waves.  This  time-domain  method 
provides  estimates  of  sea  surface  elevation  and  wave  kinematics  throughout  the  water 
column  without  compromising  the  essential  nonlinear  character  of  the  waves.  The  key  to 
the  formulation  is  retention  of  the  fully  nonlinear  free  surface  boundary  conditions 
which  are  used  in  seeking  a  nonlinear  least-squares  solution  in  a  narrow  time  window. 
Also  included  is  the  nonlinear  interaction  between  the  waves  and  a  steady  and 
depth-uniform  current.  The  analysis  method  provides  nearly  exact  results  for  theoretical 
nonlinear  steady  waves.  Application  to  field  data  requires  overspecification  in  the 
solution  to  allow  for  instrument  error  bands.  Nevertheless,  credible  solutions  were 
obtained  using  PUV  measurement  from  two  field  locations  with  very  energetic  wave  and 
current  conditions. 

Results  from  the  locally  nonlinear  analysis  differ  substantially  from  traditional 
directional  spectra.  The  method  provides  time  series  of  sea  surface  elevation  17,  radial 
frequency  w,  directional  wave  numbers  kj  and  kj,  and  Fourier  coefficients,  Aj.  These 
local  solution  parameters  provide  predictions  of  corresponding  time  series  for  any 
kinematic  quantity  at  any  position  in  the  water  column.  A  complete  field  solution  is 
available  at  each  time  interval. 

The  locally  nonlinear  PUV  method  provides  a  distinct  improvement  for  analyzing 
practical  problems  in  which  individual  wave  characteristics  are  critical.  For  example, 
wave  loads  on  offshore  structures  depend  on  wave  height,  fluid  velocity,  fluid  accelera¬ 
tion,  and  the  phasing  between  maximum  velocity  and  acceleration.  Retaining  the 
nonlinear  characteristics  inherent  in  the  measured  data  is  critical  for  credible  load 
determination. 

The  locally  nonlinear  analysis  is  fairly  computationally  intensive,  and  it  creates 
additional  information  rather  than  condensing  the  measurements.  Therefore,  this  method 
is  probably  not  well-suited  for  routine  analysis  and  archiving  of  PUV  records  obtained 
over  a  lengthy  deployment.  However,  because  the  method  does  not  rely  on  the  entire 
time  series,  it  is  possible  to  perform  selective  analysis  of  individual  waves  or  large  wave 
groups  within  the  time  series  that  may  be  of  particular  interest.  Examples  include 
so-called  ‘freak  waves’  and  even  tsunami  waves.  Also,  the  local  nature  of  the  nonlinear 
analysis  method  may  find  application  providing  near  real-time  wave  estimates  from 
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PUV  gauges.  It  may  also  assist  in  the  design  evolution  of  PUV  gauges,  as  it  provides  a 
more  realistic  interpretation  of  the  strong  nonlinearities  associated  with  the  more 
extreme  events. 
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